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ABSTRACT 

In this thesis we study several computer implementa- 
tions of the thinning algorithm, a new method for generating 
non-homogeneous Poisson processes. The method, developed 
by Professor P.A.W. Lewis, Naval Postgraduate School, 
Monterey, California, and G.S. Shedler, IBM Research 
Laboratory, San Jose, California, is valid for Poisson 
processes with any given intensity function. The basic 
thinning algorithm is modified to exploit several refine- 
ments which reduce computer execution time by approximately 
one-third. The basic and modified thinning programs are 
compared with a previous algorithm of Lewis and Shedler, 
the Poisson decomposition and gap-statistics algorithm, 

which is easily implemented for Poisson processes with 

2 

intensity functions of the form exp (aQ+a^t+a2t ). The 
thinning programs are competitive in both execution time 
and computer memory requirements. One program implementa- 
tion generates the events in a Poisson process one at a 
time; another program implements the algorithmic refinements 
which improve efficiency. 
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I . 



INTRODUCTION 



The Poisson process is a widely known and studied 
stochastic process. It is frequently used to describe 
random arrivals at some type of service facility such as 
a service station fuel pump or a bank teller's window. 

In its most common form, the "rate" of these arrivals is 
considered to be constant over time. This is the homogeneous 
Poisson process which has the familiar property that times 
between arrivals (or events) are exponentially distributed 
with mean equal to the inverse of the rate. 

The assumption of a constant rate, or homogeneity, is 
at best tenuous when applied to real world data. For exam- 
ple, the rate of arrivals at a traffic light typically 
varies from very high during rush hours to very low in the 
early morning. In addition to this cyclic time-of-day 
effect, arrival rates may exhibit longer term increases or 
decreases. Further, these effects may be superimposed upon 
shorter term effects to produce a more complex rate which 
varies with time. These processes for which arrival rates 
vary with time may often be represented by a non- homogeneous 
Poisson process, that is, a Poisson process with a time 
dependent rate of arrival. 

The generic term of Poisson process includes then both 
homogeneous and non-homogeneous Poisson processes. LEWIS 
and SHEDLER [Ref. 1] define the Poisson process generally in 
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terms of a monotone non-decreasing right-continous function 

A(t) which is bounded in any finite interval. Then the 

number of points, N(t'',t')f in any finite interval has a 

Poisson distribution with parameter y(t'',t') = A(t') — A ( t ' ' ) 

Thus, for example in (0,t'], with t' _> 0, P{N(t'',t') = n} 
n _y n 

= P{N t ,=n} = y Q e /n! , where y Q = y(0,t') = A(t') - A(0). 

The right derivative A(t) of A(t) will be assumed to 
exist and is called the rate function or intensity function 
of the process. A ( t ) is called the integrated rate function 
and has the interpretation that for t _> 0 , A(t) - A(0) = E[N ] 
For the homogeneous Poisson process A(t) is a constant, 
e.g. A, and thus the integrated rate function is simply the 
product of A and t, i.e. the expected value of = N(0,t). 

While simulation of homogeneous Poisson processes is 
relatively straightforward, the non-homogeneous Poisson 
process is more problematical. Times between events are 
not exponential in the general case and simulation has 
typically been tailored to specific classes of intensity 
functions. LEWIS and SHEDLER [Ref. 1] list three general 
methods for simulating non-homogeneous Poisson processes 
and one method for a special rate function. The general 
methods include the time scale transformation method and 
the conditioning and order-statistics method. The special 
method is the gap-statistics method, a method which is 
particular to the degree-one exponential polynomial intensity 
fucntion, i.e. those of the form A (t) = exp(b^ + b^t) . 
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Implementation of the general methods on a computer may 
pose special problems. Often the inverse of the integrated 
rate function is not explicit and must be computed numerically. 
Other problems in implementation generally result in lower 
efficiency, as measured by execution time or computer storage 
requirements or both. 

One class of intensity functions which is of general 

interest is the degree-two exponential polynomial family. 

That is, those with intensity function of the form 

2 

A(t) = exp(a^ + a^t + a^t ) . This family of functions has 
the property of being positive for all values of t, a 
necessary condition for an intensity function. Additionally, 
by varying the magnitude and sign of the coefficients, the 
exponential polynomial of degree two can be made to be mono- 
tone increasing or decreasing over time, as well as increasing 
and then decreasing, or vice versa. Use of this type of 
intensity function also leads to statistical procedures 
which are relatively simple. 

LEWIS and SHEDLER [Ref. 2] proposed a new method of 
generating the non-homogeneous Poisson process with degree- 
two exponential polynomial intensity function. It involves 
decomposition of the degree-two exponential polynomial 
intensity function, A(t), into two functions, a degree-one 
exponential polynomial function, X (t) , and a difference 
function, A_(c) = A (t) - A T (t) . This procedure allows the 
points in the degree-one exponential polynomial event stream 
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to be generated using the gap-statistics method, which is 
highly efficient when implemented on a computer. The 
remaining points with intensity function A (t) are generated 
by other methods and then merged with the other events. 

PATROW [Ref. 3] implemented two algorithms, the time 
scale transformation algorithm and the Poisson decomposi- 
tion and gap-statistics technique, and compared them for 
computational speed and computer memory requirements . His 
results indicated that the Poisson decomposition and gap- 
statistics technique was from two to seven times faster than 
the time scale transformation algorithm, although the former 
required about thirty percent more computer memory. 

PATROW' s work [Ref. 3] is also an excellent self- 
contained reference on Poisson processes, bringing many 
references together under one cover. 

A recent result of LEWIS and SHEDLER [Ref. 1] develops 
a new method for generation of points in a non-homogeneous 
Poisson process. This method, called "thinning", is similar 
to the general conditioning-acceptance-rejection method 
but has subtle differences which are computationally signi- 
ficant. The thinning method is straightforward in both an 
analytical and a computational sense, and is valid for any 
type of intensity function. The thinning theorem is 
presented in Section II. 

This thesis is, in a sense, a sequel to PATROW 's work 
[Ref. 3]. Its purpose is to implement the thinning algorithm 
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in computer program form and to compare it to the Poisson 
decomposition and gap-statistics algorithm implemented by 
PATROW [Ref. 3]. The latter implementation was designed 
for a specific subset of intensity functions, degree-two 
exponential polynomials. Since the Poisson decomposition 
and gap-statistics method outperformed a general case 
algorithm (time scale transformation) by a significant 
margin, comparing the thinning method to the Poisson 
decomposition and gap-statistics method should give a 
reasonable indication of the thinning algorithm's performance 
in generating non-homogeneous Poisson processes with other 
than degree-two exponential polynomial intensity functions. 

Section III lists the two algorithms considered, as well 
as a special application of the thinning process which will 
be of interest to those involved in event-step simulation. 
Section IV describes the methodology used in comparing the 
algorithms while Section V deals with aspects of the thinning 
procedure which may be exploited to enhance its overall 
effectiveness in a variety of situations. Finally, Section 
VI presents the results and conclusions of the comparisons 
of the algorithms. Appendices A and B contain secondary 
results and computer program listings following the 
appendices . 



